Identifying the “scaring away” pattern of wild boars from GPS trajectories in the Swiss Plateau
Ricarda Ferrari, Michael Niederberger
show data
head(wildschwein_BE)## # A tibble: 6 x 8
## TierID TierName CollarID DatetimeUTC E N day
## <int> <chr> <int> <dttm> <dbl> <dbl> <chr>
## 1 1 Ueli 12272 2014-05-28 21:01:14 2570390. 1204820. Tag
## 2 1 Ueli 12272 2014-05-28 21:15:18 2570389. 1204826. Abenddaemmerung
## 3 1 Ueli 12272 2014-05-28 21:30:13 2570391. 1204821. Abenddaemmerung
## 4 1 Ueli 12272 2014-05-28 21:45:11 2570388. 1204826. Abenddaemmerung
## 5 1 Ueli 12272 2014-05-28 22:00:33 2570388. 1204819. 1Nachtviertel
## 6 1 Ueli 12272 2014-05-28 22:15:16 2570384. 1204828. 1Nachtviertel
## # … with 1 more variable: moonilumination <dbl>head(wildschwein_metadata)## TierID TierName CollarID Sex Gewicht Study_area
## 1 1 Ueli 12272 m 79.5 Bern
## 2 1 Ueli 12844 m 91.0 Bern
## 3 2 Sabine 12275 f 62.0 Bern
## 4 5 Nicole 12273 f 50.0 Bern
## 5 10 Caroline 13570 f 68.0 Bern
## 6 10 Caroline 13969 f 58.0 Bernhead(wildschwein_overlap_temp)## # A tibble: 6 x 4
## TierID TierName CollarID Groups
## <int> <chr> <int> <dbl>
## 1 1 Ueli 12272 1
## 2 2 Sabine 12275 2
## 3 5 Nicole 12273 2
## 4 10 Caroline 13969 2
## 5 11 Isabelle 12274 2
## 6 16 Rosa 13972 2head(schreck_agenda)## # A tibble: 6 x 9
## id datum_on datum_off modus lautstaerke intervall
## <chr> <dttm> <dttm> <chr> <dbl> <dbl>
## 1 WSS_201… 2014-04-01 00:00:00 2014-06-20 00:00:00 standa… 100 15
## 2 WSS_201… 2014-07-23 00:00:00 2014-09-19 00:00:00 standa… 100 15
## 3 WSS_201… 2014-04-26 00:00:00 2014-08-08 00:00:00 standa… 50 4
## 4 WSS_201… 2014-04-26 00:00:00 2014-08-08 00:00:00 standa… 100 15
## 5 WSS_201… 2014-09-19 00:00:00 2014-10-18 00:00:00 standa… 50 20
## 6 WSS_201… 2014-05-01 00:00:00 2014-10-28 00:00:00 standa… 33 15
## # … with 3 more variables: ausrichtung_min <int>, ausrichtung_max <int>,
## # phase <dbl>head(schreck_locations)## # A tibble: 6 x 9
## id region flurname kultur installationsho… zaun jagddruck lat lon
## <chr> <chr> <chr> <chr> <dbl> <chr> <chr> <dbl> <dbl>
## 1 WSS_2… elfing… rüti weizen&h… 1.8 nein mittel 47.5 8.10
## 2 WSS_2… elfing… steiach… weizen 1.95 nein mittel 47.5 8.12
## 3 WSS_2… elfing… schlott… weizen&s… 1.8 nein mittel 47.5 8.11
## 4 WSS_2… fanel tannenh… kartoffe… 1.8 ja gering 47.0 7.06
## 5 WSS_2… fanel tannenh… karotten… 1.8 nein gering 47.0 7.06
## 6 WSS_2… fanel fanelac… kartoffe… 1.8 nein gering 47.0 7.04Schreck prepare
# Convert CRS
schreck_locations_ch <- schreck_locations %>% st_as_sf(coords = c("lon", "lat"), crs = CRS("+init=epsg:4326"), remove = FALSE)
schreck_locations_ch <- schreck_locations_ch %>% st_transform(crs = 2056)
# Filter scare device regions
schreck_locations_ch <- schreck_locations_ch %>% filter(lat < 47.2 & lon < 7.5)
coordsne <- unlist(st_geometry(schreck_locations_ch)) %>% matrix(ncol=2,byrow=TRUE) %>% as_tibble() %>% setNames(c("N","E"))
schreck_locations_ch$N <- coordsne$E
schreck_locations_ch$E <- coordsne$N
#join
schreck_locations_ch <- schreck_locations_ch %>% left_join(schreck_agenda, by=c("id"="id"))
schreck_locations_ch$wid <- c(1:25)
schreck_locations_ch <- schreck_locations_ch %>% mutate(wid=as.character(wid))Data exploratio. Get common data, show examples
# Get common samples
sabi <- wildschwein_BE %>% filter(TierName=="Sabine")
# Filter night data
sabi <- sabi %>% filter(day != "Tag" & !is.na(day)) # only night gps
s1 <- schreck_locations_ch[5,]
sabi <- sabi %>% filter(DatetimeUTC > first(s1$datum_on) & DatetimeUTC < first(s1$datum_off))
sabi <- sabi %>% mutate(dist = sqrt((first(s1$N)-N)^2+(first(s1$E)-E)^2))
sabi <- sabi %>% filter(dist < 400)
sabi <- sabi %>% mutate(trip = ifelse(hour(DatetimeUTC) > 16, day(DatetimeUTC)+1, day(DatetimeUTC)))
ggplot() +
geom_path(data=sabi %>% filter(trip<11), aes(x=E, y=N, color=factor(trip))) +
geom_point(data = s1, aes(x=E, y=N), colour="black", size=2) 
Note that the echo = FALSE parameter was added to the code chunk to prevent printing of the R code that generated the plot.
Adding closest active schreck (+ distance) to each observation
### data
w<-wildschwein_BE %>% mutate(date = as.Date(DatetimeUTC),
time= format(DatetimeUTC, format = "%H:%M:%S"),
hour= as.integer(format(DatetimeUTC, format = "%H"))+
as.integer(format(DatetimeUTC, format = "%M"))/60)
sl<-schreck_locations
###
s<-data.frame(schreck_locations_ch) ## I had to transform it to a data frame, otherwise something was weird about the coordinates (geometry)
###
s1<-s[!duplicated(s$id),c("id","N","E")] # list with all Schrecks and their location
# Set id to w
w <- w %>% mutate(id = row_number())
## reduce data to certain time frame
w1<-w %>%
filter (day!="Tag"&!is.na(day) &
((DatetimeUTC > as.Date("2014-05-01") & DatetimeUTC < as.Date("2014-07-04")) |
(DatetimeUTC > as.Date("2015-05-20") & DatetimeUTC < as.Date("2015-07-01")) |
(DatetimeUTC > as.Date("2016-04-04") & DatetimeUTC < as.Date("2016-10-04")) |
(DatetimeUTC > as.Date("2017-04-26") & DatetimeUTC < as.Date("2017-11-18"))))
nrow(w1)## [1] 38039#for(j in 1:nrow(w1)){
### look only at Schrecks that were active on that day
# s_on<-s[s$datum_on < w1[j,]$DatetimeUTC & s$datum_off > w1[j,]$DatetimeUTC,]
#if(nrow(s_on)==0)
# {w1[j,"closest_schreck"]<-"no_Schreck_on"}
#else{
#for(i in 1:nrow(s_on)){ ## calculate difference between current observation (j) and each schreck location
#s_on[i,"diff"]<-sqrt((w1[j,"N"]-s_on[i,"N"])^2+(w1[j,"E"]-s_on[i,"E"])^2)} ## add difference of current observation to location into file
#### look at the distance of the closest Schreck, only use it if less than 400m
#if(min(s_on$diff)>400){
#w1[j,"closest_schreck"]<-"no_Schreck_witin_400m"}
#else{
#w1[j,"closest_schreck"]<-s_on[s_on$diff==min(s_on$diff),"id"] ## add closest schreck to each wild boar location
#w1[j,"distance_to_closest_schreck"]<-s_on[s_on$diff==min(s_on$diff),"diff"]
#}}}
# Save data.frame to spare time
#write.csv(w1, "wildboar_loop.csv")
w1 <- read_delim("wildboar_loop.csv",",")
w1 <- w1 %>% mutate(hour= as.integer(format(DatetimeUTC, format = "%H")),
tripdate= ifelse(hour < 12, as.Date(DatetimeUTC)-1, as.Date(DatetimeUTC)),
tripdate2 = as.Date(tripdate, origin="1970-01-01"))
# Merge close wildschweinschreck gps data to origin wildboar
w <- w %>% left_join(w1 %>% dplyr::select(id, closest_schreck, distance_to_closest_schreck), by="id")
w <- w %>% mutate(closest_schreck=
ifelse(is.na(closest_schreck), "no_Schreck_witin_400m", closest_schreck))
#add tripdate
w <- w %>% mutate(hour= as.integer(format(DatetimeUTC, format = "%H")),
tripdate= ifelse(hour < 12, as.Date(DatetimeUTC)-1, as.Date(DatetimeUTC)),
tripdate2 = as.Date(tripdate, origin="1970-01-01"))
w %>% dplyr::select(DatetimeUTC, tripdate2)## # A tibble: 327,255 x 2
## DatetimeUTC tripdate2
## <dttm> <date>
## 1 2014-05-28 21:01:14 2014-05-28
## 2 2014-05-28 21:15:18 2014-05-28
## 3 2014-05-28 21:30:13 2014-05-28
## 4 2014-05-28 21:45:11 2014-05-28
## 5 2014-05-28 22:00:33 2014-05-28
## 6 2014-05-28 22:15:16 2014-05-28
## 7 2014-05-28 22:30:14 2014-05-28
## 8 2014-05-28 22:45:09 2014-05-28
## 9 2014-05-28 23:00:12 2014-05-28
## 10 2014-05-28 23:15:08 2014-05-28
## # … with 327,245 more rowsSeparation of the wildboar trips
# split time and day in separate columns:
w1$Date <- as.Date(w1$DatetimeUTC)
## order data frame by animal and time
w1<-w1[order(w1$TierName) & order(w1$DatetimeUTC),]
##### Create trips: with shorest distance to schreck (for each tripdate) ########################################
# split time and day in separate columns:
w1$Date <- as.Date(w1$DatetimeUTC)
## order data frame by animal and time
w1<-w1[order(w1$TierName) & order(w1$DatetimeUTC),]
## empty column for trip id
w$tripID_dist<-"NA"
w$start_dist<-"no"
w$isNearestPoint <- "no"
animals<-unique(w1$TierName)
for(j in 1:length(animals)){
days<-unique(w1[w1$TierName==animals[j] & !is.na(w1$TierName),]$tripdate2)
if(length(days)!=0){
for(i in 1:length(days)){
n<-w1[w1$TierName==animals[j] & w1$tripdate2==days[i],]
if(nrow(n[!is.na(n$distance_to_closest_schreck),])==0){nmin<-NA} else{
nmin<-min(n$distance_to_closest_schreck,na.rm=T)}
if(!is.na(nmin)){ ## only continue if nmin is not NA
if(nmin<=400){ ## only continue if minimal distance to schreck is less than 400m
n<-n[n$distance_to_closest_schreck==nmin & !is.na(n$distance_to_closest_schreck),]
w[w$TierName==animals[j] & w$DatetimeUTC == n$DatetimeUTC, "isNearestPoint"] <- "yes"
w[w$TierName==animals[j] & (w$DatetimeUTC >= n$DatetimeUTC-4*60*60) & (w$DatetimeUTC <= n$DatetimeUTC+4*60*60),"tripID_dist"]<-paste(animals[j],i, sep="_")
w[w$TierName==animals[j] & w$DatetimeUTC==n$DatetimeUTC,"start_dist"]<-"yes"
}}}}}
# Only trips in the night
w <- w %>% group_by(TierID, tripdate2) %>% mutate(isTrip= max(tripID_dist) != 'NA',
tripIDnight = ifelse((hour < 9.1 | hour >= 17.9)&isTrip,
max(tripID_dist), NA)) %>% ungroup()
# length(unique(w$tripIDnight)); table(w[!is.na(w$tripIDnight),]$isNearestPoint)
#List with closest_distance to schreck for each trip
trip_list <- w %>%
filter (isNearestPoint=="yes" & !is.na (tripIDnight))
nrow(trip_list)## [1] 113head(data.frame(trip_list))## TierID TierName CollarID DatetimeUTC E N day
## 1 1 Ueli 12272 2014-05-29 01:00:29 2570591 1205129 2Nachtviertel
## 2 1 Ueli 12272 2014-05-30 01:00:15 2570746 1204849 2Nachtviertel
## 3 1 Ueli 12272 2014-05-31 02:45:12 2570808 1204667 3Nachtviertel
## 4 1 Ueli 12272 2014-06-01 02:30:08 2570881 1204652 3Nachtviertel
## 5 1 Ueli 12272 2014-06-02 02:30:09 2570976 1204687 3Nachtviertel
## 6 1 Ueli 12272 2014-06-03 01:30:11 2570976 1204687 3Nachtviertel
## moonilumination date time hour id closest_schreck
## 1 0.002343097 2014-05-29 01:00:29 1 17 WSS_2014_04
## 2 0.023774932 2014-05-30 01:00:15 1 69 WSS_2014_05
## 3 0.066684515 2014-05-31 02:45:12 2 128 WSS_2014_05
## 4 0.129136625 2014-06-01 02:30:08 2 179 WSS_2014_05
## 5 0.208314672 2014-06-02 02:30:09 2 231 WSS_2014_05
## 6 0.300647726 2014-06-03 01:30:11 1 279 WSS_2014_05
## distance_to_closest_schreck tripdate tripdate2 tripID_dist start_dist
## 1 350.7307 16218 2014-05-28 Ueli_1 yes
## 2 361.1836 16219 2014-05-29 Ueli_2 yes
## 3 364.5964 16220 2014-05-30 Ueli_3 yes
## 4 318.3504 16221 2014-05-31 Ueli_4 yes
## 5 230.2212 16222 2014-06-01 Ueli_5 yes
## 6 230.8985 16223 2014-06-02 Ueli_6 yes
## isNearestPoint isTrip tripIDnight
## 1 yes TRUE Ueli_1
## 2 yes TRUE Ueli_2
## 3 yes TRUE Ueli_3
## 4 yes TRUE Ueli_4
## 5 yes TRUE Ueli_5
## 6 yes TRUE Ueli_6Overview data
### plot
pa<-ggplot(w,aes(x=date,y=TierName))+geom_line()
pa<-pa+ylab("animal ID")+scale_x_date(breaks = "3 month", minor_breaks = "1 month")+xlab("date")+theme_bw()
pa
#ps<-ggplot(s[s$id %in% w$closest_schreck,],aes(y=id,yend=id,x=datum_on,xend=datum_off))+geom_segment()
#psWildschweinschreck
# Determine mean speaker orientation
schreck_locations_ch <- schreck_locations_ch %>%
mutate(ausrichtung_mean = (ausrichtung_max - ausrichtung_min) / 2 +ausrichtung_min,
ausrichtung_dir = ifelse(ausrichtung_mean >= 45 & ausrichtung_mean < 135, "E",
ifelse(ausrichtung_mean > 135 & ausrichtung_mean < 225, "S",
ifelse(ausrichtung_mean > 225 & ausrichtung_mean < 315, "W", "N"))))
# Show example of orientation of wildboar scare device with loudness and direction
s1 <- schreck_locations_ch[5,]
wildboar_closeup <- wildschwein_BE %>%
filter(TierName=="Sabine" & DatetimeUTC > as_datetime('2015-06-08 18:30:43') & DatetimeUTC < as_datetime('2015-06-09 10:30:43')) %>%
mutate(dist = sqrt((first(s1$N)-N)^2+(first(s1$E)-E)^2)) %>% mutate(triptime0 = as.numeric((DatetimeUTC - min(DatetimeUTC))) / 60) %>% arrange(.,triptime0)
schreck_orientation <- s1 %>% mutate(length=lautstaerke*lautstaerke/100)
rad2deg <- function(rad) {(rad * 180) / (pi)}
deg2rad <- function(deg) {(deg * pi) / (180)}
schreck_orientation$asurichtung_mean <- (schreck_orientation$ausrichtung_max-
schreck_orientation$ausrichtung_min) / 2 +
schreck_orientation$ausrichtung_min
schreck_orientation$ausrichung_meanE = schreck_orientation$E[1] + schreck_orientation$length[1] * 2*
cos(deg2rad(360-270-(schreck_orientation$asurichtung_mean[1])))
schreck_orientation$ausrichung_meanN = schreck_orientation$N[1] + schreck_orientation$length[1] *2*
sin(deg2rad(360-270-(schreck_orientation$asurichtung_mean[1])))
schreck_orientation$ausrichung_minE = schreck_orientation$E[1] + schreck_orientation$length[1] *
cos(deg2rad(360-270-schreck_orientation$ausrichtung_min[1]))
schreck_orientation$ausrichung_minN = schreck_orientation$N[1] + schreck_orientation$length[1] *
sin(deg2rad(360-270-schreck_orientation$ausrichtung_min[1]))
schreck_orientation$ausrichung_maxE = schreck_orientation$E[1] + schreck_orientation$length[1] *
cos(deg2rad(360-270-schreck_orientation$ausrichtung_max[1]))
schreck_orientation$ausrichung_maxN = schreck_orientation$N[1] + schreck_orientation$length[1] *
sin(deg2rad(360-270-schreck_orientation$ausrichtung_max[1]))
x_coord <- c(schreck_orientation$E[1], schreck_orientation$ausrichung_minE[1],
schreck_orientation$ausrichung_meanE[1],
schreck_orientation$ausrichung_maxE[1], schreck_orientation$E[1])
y_coord <- c(schreck_orientation$N[1], schreck_orientation$ausrichung_minN[1],
schreck_orientation$ausrichung_meanN[1],
schreck_orientation$ausrichung_maxN[1], schreck_orientation$N[1])
poly <- st_polygon(list(matrix(c(x_coord, y_coord),ncol=2, byrow=FALSE)))
pos <- data.frame(x=x_coord, y=y_coord, id=c(1,2,3,4,5))
# Show example of scare device orientation and loudness
ggplot() +
geom_path(data=wildboar_closeup, aes(x=E, y=N, color=triptime0)) +
geom_point(data = s1, aes(x=E, y=N), colour="black", size=2) +
geom_point(data = schreck_orientation, aes(x=ausrichung_meanE, y=ausrichung_meanN),
colour="brown", size=2) +
geom_polygon(data=pos, aes(x=x, y = y), fill="orange", alpha=0.4)
Wilddboar scare effect properties
# Get all trips
trips <- w %>% filter(!is.na(tripIDnight)) %>% group_by(tripIDnight)
trips ## # A tibble: 6,661 x 21
## # Groups: tripIDnight [113]
## TierID TierName CollarID DatetimeUTC E N day
## <int> <chr> <int> <dttm> <dbl> <dbl> <chr>
## 1 1 Ueli 12272 2014-05-28 21:01:14 2570390. 1204820. Tag
## 2 1 Ueli 12272 2014-05-28 21:15:18 2570389. 1204826. Abenddaemmeru…
## 3 1 Ueli 12272 2014-05-28 21:30:13 2570391. 1204821. Abenddaemmeru…
## 4 1 Ueli 12272 2014-05-28 21:45:11 2570388. 1204826. Abenddaemmeru…
## 5 1 Ueli 12272 2014-05-28 22:00:33 2570388. 1204819. 1Nachtviertel
## 6 1 Ueli 12272 2014-05-28 22:15:16 2570384. 1204828. 1Nachtviertel
## 7 1 Ueli 12272 2014-05-28 22:30:14 2570393. 1204824. 1Nachtviertel
## 8 1 Ueli 12272 2014-05-28 22:45:09 2570585. 1205044. 1Nachtviertel
## 9 1 Ueli 12272 2014-05-28 23:00:12 2570576. 1205044. 1Nachtviertel
## 10 1 Ueli 12272 2014-05-28 23:15:08 2570566. 1205047. 1Nachtviertel
## # … with 6,651 more rows, and 14 more variables: moonilumination <dbl>,
## # date <date>, time <chr>, hour <int>, id <dbl>, closest_schreck <chr>,
## # distance_to_closest_schreck <dbl>, tripdate <dbl>, tripdate2 <date>,
## # tripID_dist <chr>, start_dist <chr>, isNearestPoint <chr>, isTrip <lgl>,
## # tripIDnight <chr>wildboar_trip_scared <- data.frame(matrix(ncol = 22, nrow = 0))
x <- c("id", "TierName","triptime0", "approachingRate", "approachingRateRelative",
"approachingRateAbsolute", "speed", "sinousity", "linedist", "acceleration", "speedDiff3",
"dist","E", "N","x", "y", "directionRelative",
"closest_schreck", "tripIDnight", "DatetimeUTC", "hour", "day")
colnames(wildboar_trip_scared) <- x
# Amount of trips
unique(trips$tripIDnight)## [1] "Ueli_1" "Ueli_2" "Ueli_3" "Ueli_4" "Ueli_5" "Ueli_6"
## [7] "Ueli_7" "Ueli_9" "Ueli_10" "Ueli_12" "Ueli_14" "Ueli_15"
## [13] "Ueli_110" "Ueli_114" "Ueli_119" "Ueli_124" "Ueli_127" "Ueli_133"
## [19] "Ueli_137" "Sabine_2" "Sabine_18" "Sabine_20" "Sabine_21" "Sabine_22"
## [25] "Sabine_23" "Sabine_24" "Sabine_25" "Sabine_26" "Sabine_27" "Sabine_28"
## [31] "Sabine_29" "Sabine_30" "Sabine_31" "Sabine_32" "Sabine_33" "Sabine_34"
## [37] "Sabine_35" "Sabine_36" "Sabine_37" "Sabine_38" "Sabine_39" "Sabine_40"
## [43] "Sabine_41" "Sabine_42" "Sabine_43" "Ruth_6" "Ruth_7" "Ruth_8"
## [49] "Ruth_9" "Ruth_10" "Ruth_12" "Ruth_14" "Ruth_15" "Ruth_16"
## [55] "Ruth_17" "Ruth_18" "Ruth_19" "Ruth_20" "Ruth_21" "Ruth_22"
## [61] "Ruth_23" "Ruth_24" "Ruth_25" "Ruth_26" "Ruth_27" "Ruth_28"
## [67] "Ruth_29" "Ruth_30" "Ruth_31" "Ruth_32" "Ruth_33" "Ruth_34"
## [73] "Ruth_35" "Ruth_36" "Ruth_37" "Ruth_38" "Ruth_39" "Ruth_40"
## [79] "Ruth_41" "Ruth_42" "Ruth_43" "Olga_1" "Olga_2" "Olga_4"
## [85] "Olga_6" "Olga_11" "Olga_13" "Olga_14" "Olga_15" "Olga_16"
## [91] "Olga_18" "Olga_19" "Olga_20" "Olga_21" "Olga_22" "Olga_23"
## [97] "Olga_24" "Olga_25" "Olga_26" "Olga_27" "Olga_28" "Olga_29"
## [103] "Olga_31" "Olga_33" "Olga_34" "Olga_35" "Olga_36" "Olga_37"
## [109] "Olga_38" "Olga_40" "Olga_41" "Olga_42" "Olga_43"for (u in 1:length(unique(trips$tripIDnight)))
{
wildboar_trip <- trips %>% filter(tripIDnight==unique(trips$tripIDnight)[u])
# Only take full trips
if (length(wildboar_trip$TierID) > 61)
{
# Get wildboard schreck associated to wild boar
schreck <- schreck_locations_ch %>% filter(id == max(wildboar_trip$closest_schreck))
# Calculate distance and time differences
wildboar_trip <- wildboar_trip %>%
mutate(dist = sqrt((first(schreck$N)-N)^2+(first(schreck$E)-E)^2),
triptime0 = as.numeric((DatetimeUTC - min(DatetimeUTC))) / 60,
triptimeDiff = (triptime0-lag(triptime0))*60,
distanceAbsolute = ((E- lag(E))^2 + (N-lag(N))^2)^0.5) %>% arrange(.,triptime0)
# Calculate approaching rate
wildboar_trip <- wildboar_trip %>%
mutate(approachingRate = lag(dist)-dist,
approachingRateAbsolute = approachingRate / (triptimeDiff),
approachingRateRelative = approachingRate / (distanceAbsolute),
approachingRateRelative = ifelse(is.na(approachingRateRelative), 0, approachingRateRelative))
# Calculate speed & acceleration
wildboar_trip <- wildboar_trip %>%
mutate(speed = round(distanceAbsolute / triptimeDiff, 4),
acceleration = round((lead(speed)-speed) * 60 / lead(triptimeDiff), 4),
speedDiff3 = round(((speed+lead(speed)+lead(speed, 2))/3 -
(lag(speed)+lag(speed, 2)+lag(speed, 3)/3)), 4))
# Relative coordinates (trajectory)
wildboar_trip <- wildboar_trip %>% mutate(x=E-first(schreck$E), y=N-first(schreck$N))
coords <- data.frame(x = wildboar_trip$x,
y = wildboar_trip$y,
times = wildboar_trip$triptime0)
# Set all sinousitiy to 0
wildboar_trip$sinousity = replicate(length(coords$x), NA)
for (i in 1:(length(coords$x)-3)) {
# Create a trajectory from the coordinates
trj <- TrajFromCoords(coords[(i):(i+3),])
# Rescale stepsize
trj_re <- TrajRediscretize(trj, 1)
# Calculate sinousity
wildboar_trip$sinousity[i] = TrajSinuosity(trj_re, compass.direction = TRUE) %>%
round(digits = 4)
}
# Get proximity of movement
wildboar_trip <- wildboar_trip %>%
mutate(linedist = ((x^2+y^2)^0.5 + ((x-lag(x))^2+(y-lag(y))^2)^0.5)/2)
# Get relative direction to scare device
wildboar_trip <- wildboar_trip %>% mutate(directionRelative=
ifelse(abs(x) > abs(y),
ifelse(x>0, "E", "W"),
ifelse(y>0, "N", "S")))
# Select
wildboar_trip_scared <- wildboar_trip_scared %>%
rbind(wildboar_trip %>%
dplyr::select(id, TierName, triptime0, approachingRate, approachingRateRelative,
approachingRateAbsolute, speed, sinousity, linedist, acceleration,
speedDiff3, dist, E, N, x, y, directionRelative, closest_schreck, tripIDnight,
DatetimeUTC, hour, day))
}
}
head(wildboar_trip_scared)## # A tibble: 6 x 22
## # Groups: tripIDnight [1]
## id TierName triptime0 approachingRate approachingRateRe… approachingRateAb…
## <dbl> <chr> <dbl> <dbl> <dbl> <dbl>
## 1 15062 Ueli 0 NA 0 NA
## 2 15063 Ueli 15.1 3.68 0.717 0.00407
## 3 15064 Ueli 30.0 3.38 0.488 0.00376
## 4 15065 Ueli 45.1 0.606 0.111 0.000670
## 5 15066 Ueli 60.6 -6.72 -0.901 -0.00726
## 6 15067 Ueli 75.1 9.57 0.869 0.0110
## # … with 16 more variables: speed <dbl>, sinousity <dbl>, linedist <dbl>,
## # acceleration <dbl>, speedDiff3 <dbl>, dist <dbl>, E <dbl>, N <dbl>,
## # x <dbl>, y <dbl>, directionRelative <chr>, closest_schreck <chr>,
## # tripIDnight <chr>, DatetimeUTC <dttm>, hour <int>, day <chr># Scaled values, calculate scare effect and set "scaring away" points
wildboar_trip_scared <- wildboar_trip_scared %>% group_by(tripIDnight) %>%
mutate(
approachingRateRelativeS = (approachingRateRelative-min(approachingRateRelative, na.rm = TRUE)) /
(max(approachingRateRelative, na.rm = TRUE)- min(approachingRateRelative, na.rm = TRUE)),
accelerationS = (acceleration-min(acceleration, na.rm = TRUE))/
(max(acceleration, na.rm = TRUE) - min(acceleration, na.rm = TRUE)),
sinousityS = (sinousity-min(sinousity, na.rm = TRUE))/
(max(sinousity, na.rm = TRUE)- min(sinousity, na.rm = TRUE)),
distanceS = (dist-min(dist, na.rm = TRUE)) / (max(dist, na.rm = TRUE) - min(dist, na.rm = TRUE)),
speedDiff3S = (speedDiff3-min(speedDiff3, na.rm=TRUE)) /
(max(speedDiff3, na.rm=TRUE)-min(speedDiff3, na.rm = TRUE)),
approachingRateRelativeSlead = lead(approachingRateRelativeS),
scareEffect = accelerationS - sinousityS - distanceS - lead(approachingRateRelativeS) + speedDiff3S,
scareEffectShow = ifelse(scareEffect < 0, -0.1, scareEffect),
scary= scareEffect > 1.4,
scary= ifelse(is.na(scary), FALSE, scary),
scary= ifelse(scary&lag(scary), FALSE, scary)) %>%
ungroup()
# Get scared points
wildboar_scared <- wildboar_trip_scared %>% filter(scary)
# Get scared trips
wildboar_trip_scared_True <- wildboar_trip_scared %>%
group_by(tripIDnight) %>%
mutate(tripScared = sum(scary, na.rm = TRUE)) %>%
filter(tripScared > 0)
# Join trips with scared wild boars
w <- w %>% left_join(wildboar_trip_scared_True %>% ungroup() %>%
dplyr::select(id, scareEffect, triptime0, accelerationS,
approachingRateRelativeSlead, speedDiff3S, sinousityS,
distanceS, x, y, scary), by="id")Plot all found shrecks
# Amount of scared trips
unique(wildboar_trip_scared$tripIDnight)## [1] "Ueli_110" "Ueli_114" "Ueli_119" "Ueli_124" "Ueli_127" "Ueli_133"
## [7] "Ueli_137" "Sabine_2" "Sabine_18" "Sabine_20" "Sabine_21" "Sabine_22"
## [13] "Sabine_23" "Sabine_24" "Sabine_25" "Sabine_26" "Sabine_27" "Sabine_28"
## [19] "Sabine_29" "Sabine_30" "Sabine_31" "Sabine_32" "Sabine_33" "Sabine_34"
## [25] "Sabine_35" "Sabine_36" "Sabine_37" "Sabine_38" "Sabine_39" "Sabine_40"
## [31] "Sabine_41" "Sabine_42" "Sabine_43" "Ruth_12" "Ruth_17" "Ruth_40"
## [37] "Ruth_42" "Ruth_43" "Olga_1" "Olga_2" "Olga_4" "Olga_6"
## [43] "Olga_11" "Olga_13" "Olga_14" "Olga_15" "Olga_16" "Olga_18"
## [49] "Olga_19" "Olga_20" "Olga_21" "Olga_22" "Olga_23" "Olga_24"
## [55] "Olga_25" "Olga_26" "Olga_27" "Olga_28" "Olga_29" "Olga_31"
## [61] "Olga_33" "Olga_34" "Olga_35" "Olga_36" "Olga_37" "Olga_38"
## [67] "Olga_40" "Olga_41" "Olga_42" "Olga_43"unique(wildboar_trip_scared_True$tripIDnight)## [1] "Sabine_21" "Sabine_23" "Sabine_24" "Sabine_28" "Sabine_31" "Sabine_35"
## [7] "Sabine_36" "Sabine_39" "Sabine_40" "Sabine_43" "Ruth_40" "Ruth_43"
## [13] "Olga_4" "Olga_16" "Olga_19" "Olga_24" "Olga_27" "Olga_28"
## [19] "Olga_31" "Olga_33" "Olga_34" "Olga_35" "Olga_36" "Olga_43"length(unique(wildboar_trip_scared_True$tripIDnight))## [1] 24for (u in 1:length(unique(wildboar_trip_scared_True$tripIDnight)))
{
test <- wildboar_trip_scared_True %>% filter(tripIDnight==unique(wildboar_trip_scared_True$tripIDnight)[u])
print(first(test$DatetimeUTC))
print(first(test$tripIDnight))
# Variable values
print(ggplot(data=test) +
geom_bar(stat="identity", aes(x=triptime0/60, y=scareEffectShow),fill = "grey", size=1.3)+
geom_line( aes(x=triptime0/60, y=distanceS, alpha=0.8), color="blue", size=1.3) +
geom_line( aes(x=triptime0/60, y=approachingRateRelativeSlead, alpha=0.8), color="green", size=1.3) +
geom_line( aes(x=triptime0/60, y=speedDiff3S, alpha=0.8), color="red", size=1.3) +
geom_line( aes(x=triptime0/60, y=sinousityS, alpha=0.8), color="turquoise", size=1.3) +
geom_line( aes(x=triptime0/60, y=accelerationS, alpha=0.8), color="orange", size=1.3) +
ggtitle(paste("Wildboar trip", first(test$tripIDnight), "started at", first(test$DatetimeUTC))) +
xlab("Hour passed by from 18:00 p.m.") + ylab("Scaled variables (0 to 1) and scare effect (-3 to 2)"))
testClose <- test %>%
filter(scary | lead(scary) | lead(scary,2)| lead(scary,3)| lead(scary,4) | lag(scary,1)| lag(scary,2)| lag(scary,3)| lag(scary,4))
print(first(test$DatetimeUTC))
print(first(test$tripIDnight))
print(ggplot(data=testClose) +
geom_bar(stat="identity", aes(x=(triptime0-testClose$triptime0[5]), y=scareEffectShow),
fill = "grey", size=1.3)+
geom_line(aes(x=(triptime0-testClose$triptime0[5]), y=distanceS, alpha=0.8, colour="Distance"),
color="blue", size=1.3) +
geom_line( aes(x=(triptime0-testClose$triptime0[5]),
y=approachingRateRelativeSlead, alpha=0.8), color="green", size=1.3) +
geom_line( aes(x=(triptime0-testClose$triptime0[5]), y=accelerationS, alpha=0.8), color="orange", size=1.3) +
geom_line( aes(x=(triptime0-testClose$triptime0[5]), y=speedDiff3S, alpha=0.8), color="red", size=1.3) +
geom_line( aes(x=(triptime0-testClose$triptime0[5]), y=sinousityS, alpha=0.8), color="turquoise", size=1.3) +
ggtitle(paste("Wildboar trip", first(testClose$tripIDnight), "1 hour before and after scare effect at",
testClose$DatetimeUTC[5])) +
xlab("Minutes before or after the scare effect") +
ylab("Scaled variables (0 to 1) and scare effect (-3 to 2)")+
theme(legend.position = "bottom"))
print(ggplot() +
geom_path(data= testClose, aes(x=E, y=N, color=(triptime0-testClose$triptime0[5]))) +
geom_point(data = schreck_locations_ch %>% filter(id == max(test$closest_schreck)),
aes(x=E, y=N), colour="black", size=2) +
geom_point(data = test %>% filter(scary), aes(x=E, y=N), color="red")+
ggtitle(paste("Wildboar trip", first(test$tripIDnight), "started at", first(test$DatetimeUTC))))
print(ggplot() +
geom_path(data= test, aes(x=E, y=N, color=triptime0)) +
geom_point(data = schreck_locations_ch %>% filter(id == max(test$closest_schreck)),
aes(x=E, y=N), colour="black", size=2) +
geom_point(data = test %>% filter(scary), aes(x=E, y=N), color="red")+
ggtitle(paste("Wildboar trip", first(testClose$tripIDnight), "1 hour before and after scare effect at",
testClose$DatetimeUTC[5])))
}## [1] "2015-06-08 18:01:16 UTC"
## [1] "Sabine_21"
## [1] "2015-06-08 18:01:16 UTC"
## [1] "Sabine_21"
## [1] "2015-06-10 18:00:20 UTC"
## [1] "Sabine_23"
## [1] "2015-06-10 18:00:20 UTC"
## [1] "Sabine_23"
## [1] "2015-06-11 18:00:12 UTC"
## [1] "Sabine_24"
## [1] "2015-06-11 18:00:12 UTC"
## [1] "Sabine_24"
## [1] "2015-06-15 18:00:22 UTC"
## [1] "Sabine_28"
## [1] "2015-06-15 18:00:22 UTC"
## [1] "Sabine_28"
## [1] "2015-06-18 18:00:12 UTC"
## [1] "Sabine_31"
## [1] "2015-06-18 18:00:12 UTC"
## [1] "Sabine_31"
## [1] "2015-06-22 18:00:44 UTC"
## [1] "Sabine_35"
## [1] "2015-06-22 18:00:44 UTC"
## [1] "Sabine_35"
## [1] "2015-06-23 18:00:08 UTC"
## [1] "Sabine_36"
## [1] "2015-06-23 18:00:08 UTC"
## [1] "Sabine_36"
## [1] "2015-06-26 18:00:15 UTC"
## [1] "Sabine_39"
## [1] "2015-06-26 18:00:15 UTC"
## [1] "Sabine_39"
## [1] "2015-06-27 18:00:10 UTC"
## [1] "Sabine_40"
## [1] "2015-06-27 18:00:10 UTC"
## [1] "Sabine_40"
## [1] "2015-06-30 18:01:08 UTC"
## [1] "Sabine_43"
## [1] "2015-06-30 18:01:08 UTC"
## [1] "Sabine_43"
## [1] "2015-06-27 18:02:34 UTC"
## [1] "Ruth_40"
## [1] "2015-06-27 18:02:34 UTC"
## [1] "Ruth_40"
## [1] "2015-06-30 18:00:43 UTC"
## [1] "Ruth_43"
## [1] "2015-06-30 18:00:43 UTC"
## [1] "Ruth_43"
## [1] "2015-05-22 18:00:15 UTC"
## [1] "Olga_4"
## [1] "2015-05-22 18:00:15 UTC"
## [1] "Olga_4"
## [1] "2015-06-03 18:00:16 UTC"
## [1] "Olga_16"
## [1] "2015-06-03 18:00:16 UTC"
## [1] "Olga_16"
## [1] "2015-06-06 18:00:11 UTC"
## [1] "Olga_19"
## [1] "2015-06-06 18:00:11 UTC"
## [1] "Olga_19"
## [1] "2015-06-11 18:00:20 UTC"
## [1] "Olga_24"
## [1] "2015-06-11 18:00:20 UTC"
## [1] "Olga_24"
## [1] "2015-06-14 18:00:38 UTC"
## [1] "Olga_27"
## [1] "2015-06-14 18:00:38 UTC"
## [1] "Olga_27"
## [1] "2015-06-15 18:00:21 UTC"
## [1] "Olga_28"
## [1] "2015-06-15 18:00:21 UTC"
## [1] "Olga_28"
## [1] "2015-06-18 18:00:10 UTC"
## [1] "Olga_31"
## [1] "2015-06-18 18:00:10 UTC"
## [1] "Olga_31"
## [1] "2015-06-20 18:00:13 UTC"
## [1] "Olga_33"
## [1] "2015-06-20 18:00:13 UTC"
## [1] "Olga_33"
## [1] "2015-06-21 18:00:09 UTC"
## [1] "Olga_34"
## [1] "2015-06-21 18:00:09 UTC"
## [1] "Olga_34"
## [1] "2015-06-22 18:00:10 UTC"
## [1] "Olga_35"
## [1] "2015-06-22 18:00:10 UTC"
## [1] "Olga_35"
## [1] "2015-06-23 18:00:12 UTC"
## [1] "Olga_36"
## [1] "2015-06-23 18:00:12 UTC"
## [1] "Olga_36"
## [1] "2015-06-30 18:00:37 UTC"
## [1] "Olga_43"
## [1] "2015-06-30 18:00:37 UTC"
## [1] "Olga_43"
Calculate noise exposition at scared data points and for trips.
# Get scared points and calculate sound volume
wildboar_scared_schreck <- wildboar_scared %>%
left_join(schreck_locations_ch %>% rename(schreckID = id), by=c("closest_schreck"="schreckID")) %>%
filter(datum_on < DatetimeUTC & datum_off > DatetimeUTC)
wildboar_scared_schreck <- wildboar_scared_schreck %>%
select(id, tripIDnight, TierName, DatetimeUTC, ausrichtung_dir, directionRelative, dist, lautstaerke, scareEffect, modus, day) %>%
mutate(sameDirection=ausrichtung_dir==directionRelative,
volumeDir = ifelse(sameDirection, lautstaerke,
ifelse(ausrichtung_dir=='S'&directionRelative=='N' |
ausrichtung_dir=='N'& directionRelative=='S' |
ausrichtung_dir=='W'& directionRelative=='E'|
ausrichtung_dir=='E'& directionRelative=='W',
lautstaerke-10, lautstaerke-5)),
soundDistance = volumeDir - abs(20 * log(1/(1.581*dist), 10)),
soundDistance = ifelse(day=="Tag", -20, soundDistance),
schreckProbabily = soundDistance *100/40)
ggplot(data=wildboar_scared_schreck) +
geom_point(aes(x= scareEffect, y=soundDistance, color=soundDistance>10), size=2)
# Get scared trips and calculate sound volume
wildboar_trips <- wildboar_trip_scared %>% group_by(tripIDnight) %>%
filter(scareEffect == max(scareEffect, na.rm=TRUE)) %>% ungroup() %>%
left_join(schreck_locations_ch, by=c("closest_schreck"="id")) %>%
filter(datum_on < DatetimeUTC & datum_off > DatetimeUTC)
wildboar_trips <- wildboar_trips %>%
select(id, tripIDnight, TierName, DatetimeUTC, scary, ausrichtung_dir, directionRelative,
dist, lautstaerke, scareEffect, modus, day) %>%
mutate(sameDirection=ausrichtung_dir==directionRelative,
volumeDir = ifelse(sameDirection, lautstaerke,
ifelse(ausrichtung_dir=='S' & directionRelative=='N' |
ausrichtung_dir=='N' & directionRelative=='S' |
ausrichtung_dir=='W' & directionRelative=='E'|
ausrichtung_dir=='E' & directionRelative=='W',
lautstaerke-10, lautstaerke-5)),
soundDistance = volumeDir - abs(20 * log(1/(1.581*dist), 10)),
soundDistance = ifelse(day=="Tag", -20, soundDistance),
scareSound= ifelse(soundDistance < 10, FALSE, scary))
wildboar_trips$sum <- replicate(length(wildboar_trips$TierName), 1)
# Filter scare effect trips by scare sound
wildboar_scared2 <- wildboar_scared %>%
left_join(wildboar_trips %>% select(id, soundDistance, scareSound), by=c("id"="id"))
wildboar_scared2 <- wildboar_scared2 %>% filter(scareSound)
wildboar_scared2$ScareAndSound <- replicate(length(wildboar_scared2$TierName), TRUE)
w <- w %>% left_join(wildboar_scared2 %>% select(id, soundDistance, ScareAndSound),by=c('id'='id'))
# Get scared trips and calculate sound volume, select most intense sound point per trip and highest scare effect per trip.
wildboar_trips2 <- wildboar_trip_scared %>%
left_join(schreck_locations_ch, by=c("closest_schreck"="id")) %>%
filter(datum_on < DatetimeUTC & datum_off > DatetimeUTC) %>%
select(tripIDnight, TierName, DatetimeUTC, scary, ausrichtung_dir, directionRelative,
dist, lautstaerke, scareEffect, modus, day) %>%
mutate(sameDirection=ausrichtung_dir==directionRelative,
volumeDir = ifelse(sameDirection, lautstaerke,
ifelse(ausrichtung_dir=='S' & directionRelative=='N' |
ausrichtung_dir=='N' & directionRelative=='S' |
ausrichtung_dir=='W' & directionRelative=='E'|
ausrichtung_dir=='E' & directionRelative=='W',
lautstaerke-10, lautstaerke-5)),
soundDistance = volumeDir - abs(20 * log(1/(1.581*dist), 10)),
soundDistance = ifelse(day=="Tag", -20, soundDistance),
schreckProbabily = soundDistance *100/40)
wildboar_trips2 <- wildboar_trips2 %>% group_by(tripIDnight) %>%
mutate(scareEffectMaxTrip=max(scareEffect)) %>%
filter(soundDistance == max(soundDistance, na.rm=TRUE)) %>% ungroup()
wildboar_trips2$sum <- replicate(length(wildboar_trips2$TierName), 1)
ggplot(data=wildboar_trips2) +
geom_point(aes(x= scareEffectMaxTrip, y=soundDistance, color=soundDistance>10 & scareEffectMaxTrip > 1.4), size=2)
# Show plot of trips (total, noise exposed, scared, and both)
ggplot(data=wildboar_trips) +
geom_bar(stat="identity", aes(x=TierName, y=sum),
fill = "grey", size=1.3) +
geom_bar(data = wildboar_trips2 %>% filter(soundDistance > 20), stat="identity", aes(x=TierName, y= sum),
fill = "#ff9999", size=1.3) +
geom_bar(data = wildboar_trips %>% filter(scary), stat="identity", aes(x=TierName, y=sum),
fill = "#e83c3c", size=1.3) +
geom_bar(data = wildboar_trips %>% filter(soundDistance > 10 & scary), stat="identity", aes(x=TierName, y= sum),
fill = "darkred", size=1.3) +
xlab("Wild boar animal") +
ylab("Amount of trips")
wildboar_trips %>% filter(soundDistance > 20 & scary)## # A tibble: 13 x 17
## id tripIDnight TierName DatetimeUTC scary ausrichtung_dir
## <dbl> <chr> <chr> <dttm> <lgl> <chr>
## 1 42304 Sabine_24 Sabine 2015-06-12 02:00:43 TRUE S
## 2 42689 Sabine_28 Sabine 2015-06-16 02:15:09 TRUE S
## 3 42970 Sabine_31 Sabine 2015-06-19 00:30:13 TRUE S
## 4 43745 Sabine_39 Sabine 2015-06-27 02:15:11 TRUE S
## 5 43838 Sabine_40 Sabine 2015-06-28 01:30:10 TRUE S
## 6 193706 Olga_24 Olga 2015-06-12 02:15:17 TRUE S
## 7 193979 Olga_27 Olga 2015-06-14 22:45:13 TRUE S
## 8 194085 Olga_28 Olga 2015-06-16 01:15:09 TRUE S
## 9 194373 Olga_31 Olga 2015-06-19 01:30:10 TRUE S
## 10 194558 Olga_33 Olga 2015-06-20 23:45:12 TRUE S
## 11 194664 Olga_34 Olga 2015-06-22 02:15:12 TRUE S
## 12 194768 Olga_35 Olga 2015-06-23 04:15:14 TRUE S
## 13 194845 Olga_36 Olga 2015-06-23 23:30:11 TRUE S
## # … with 11 more variables: directionRelative <chr>, dist <dbl>,
## # lautstaerke <dbl>, scareEffect <dbl>, modus <chr>, day <chr>,
## # sameDirection <lgl>, volumeDir <dbl>, soundDistance <dbl>,
## # scareSound <lgl>, sum <dbl>schreck_locations_ch## Simple feature collection with 25 features and 22 fields
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: 2549490 ymin: 1184326 xmax: 2582053 ymax: 1214438
## Projected CRS: CH1903+ / LV95
## # A tibble: 25 x 23
## id region flurname kultur installationsho… zaun jagddruck lat lon
## * <chr> <chr> <chr> <chr> <dbl> <chr> <chr> <dbl> <dbl>
## 1 WSS_2… fanel tannenhof kartoff… 1.8 ja gering 47.0 7.06
## 2 WSS_2… fanel tannenhof karotte… 1.8 nein gering 47.0 7.06
## 3 WSS_2… fanel fanelach… kartoff… 1.8 nein gering 47.0 7.04
## 4 WSS_2… fanel fanelach… kartoff… 1.8 nein gering 47.0 7.04
## 5 WSS_2… fanel tannenhof weizen 1.8 nein gering 47.0 7.06
## 6 WSS_2… fanel tannenhof weizen 1.8 nein gering 47.0 7.06
## 7 WSS_2… fanel fanelach… weizen 1.8 nein gering 47.0 7.04
## 8 WSS_2… fanel fanelach… weizen 1.8 nein gering 47.0 7.04
## 9 WSS_2… fanel gampelen… weizen&… 2 nein gering 47.0 7.05
## 10 WSS_2… fanel gampelen… weizen&… 2 nein gering 47.0 7.05
## # … with 15 more rows, and 14 more variables: geometry <POINT [m]>, N <dbl>,
## # E <dbl>, datum_on <dttm>, datum_off <dttm>, modus <chr>, lautstaerke <dbl>,
## # intervall <dbl>, ausrichtung_min <int>, ausrichtung_max <int>, phase <dbl>,
## # wid <chr>, ausrichtung_mean <dbl>, ausrichtung_dir <chr>How long does the scare-off effect last
wildboar_trip_scared$return <- replicate(length(wildboar_trip_scared$N), "NA")
wildboar_trip_scared$return <- as.character(wildboar_trip_scared$return)
wildboar_trip_scared$return_hrs <- NA
### looking after end of trip
trip_list$hour<- as.integer(format(trip_list$DatetimeUTC, format = "%H")) ## add hour variable to trip list
trip_list$date2<-trip_list$date
trip_list[trip_list$hour<12,]$date2<-trip_list[trip_list$hour<12,]$date-1
for(j in 1:nrow(trip_list)){
wildboar_trip_scared$return <- as.character(wildboar_trip_scared$return)
wildboar_trip_scared$return_hrs <- as.numeric(as.character(wildboar_trip_scared$return_hrs))
tripend<-paste(trip_list$date2[j]+1,"09:45:00 UTC",sep=" ") ## end of trip
under400 <- wildboar_trip_scared[wildboar_trip_scared$TierName==trip_list$TierName[j] & wildboar_trip_scared$DatetimeUTC>tripend & wildboar_trip_scared$closest_schreck==trip_list$closest_schreck[j],]
under400 <- under400[order(under400$DatetimeUTC),]
if(nrow(under400)==0){wildboar_trip_scared[wildboar_trip_scared$tripIDnight==trip_list$tripIDnight[j] & !is.na(wildboar_trip_scared$tripIDnight),]$return <- "never"
wildboar_trip_scared[wildboar_trip_scared$tripIDnight==trip_list$tripIDnight[j] & !is.na(wildboar_trip_scared$tripIDnight), "return_hrs"] <- NA
}
else{
wildboar_trip_scared[wildboar_trip_scared$tripIDnight==trip_list$tripIDnight[j] & !is.na(wildboar_trip_scared$tripIDnight), "return"] <- as.character(under400[1,]$DatetimeUTC)
wildboar_trip_scared$return <- as.POSIXct(wildboar_trip_scared$return, format="%Y-%m-%d %H:%M:%S",tz="UTC")
return_time <- as.numeric(as.character(difftime(under400[1,]$DatetimeUTC,trip_list[j,]$DatetimeUTC, units="hours")))
wildboar_trip_scared[wildboar_trip_scared$tripIDnight==trip_list$tripIDnight[j] & !is.na(wildboar_trip_scared$tripIDnight),"return_hrs"] <- return_time
}}
### create file with one line per trip (the observation with the highest scare away effect)
wildboar_trip_scared$highest<-"no"
triplist <- unique(wildboar_trip_scared$tripIDnight)
for(i in 1:length(triplist)){
list<-wildboar_trip_scared[wildboar_trip_scared$tripIDnight==triplist[i],"scareEffect"]
maxscare<-max(list,na.rm=T)
if(!is.na(maxscare)){
wildboar_trip_scared[wildboar_trip_scared$scareEffect==maxscare & !is.na(wildboar_trip_scared$scareEffect),"highest"]<-"yes"}}
### unique file
ww<-wildboar_trip_scared[wildboar_trip_scared$highest=="yes",]
### creating new column scared_away (yes if at least one observation is classified as scared for a given trip)
ww$scared_away <- "no"
for(i in 1:length(triplist)){
list<-ww[ww$tripIDnight==triplist[i],"scary"]
if(TRUE %in% list){ww[ww$tripIDnight==triplist[i],]$scared_away<-"yes"}}
tt<-ww
head(data.frame(tt))## id TierName triptime0 approachingRate approachingRateRelative
## 1 15095 Ueli 495.0500 175.10762 0.6963968
## 2 15462 Ueli 239.9833 4.96046 0.8264131
## 3 15959 Ueli 494.9500 -133.30475 -0.4795234
## 4 16423 Ueli 254.8333 -11.08749 -0.8571269
## 5 16711 Ueli 254.5500 -13.08219 -0.1301012
## 6 17298 Ueli 419.8333 -85.31995 -0.4489033
## approachingRateAbsolute speed sinousity linedist acceleration speedDiff3
## 1 0.195432611 0.2806 0.2239 275.01180 0.0109 0.3262
## 2 0.005536228 0.0067 0.5851 57.58407 0.0001 -0.0442
## 3 -0.148777618 0.3103 1.0479 313.01142 0.0110 0.0014
## 4 -0.012374431 0.0144 0.3875 31.11236 0.0041 -0.0594
## 5 -0.014616972 0.1124 1.1569 183.66269 0.0118 0.0341
## 6 -0.095117006 0.2119 0.1238 524.28555 0.0094 0.2170
## dist E N x y directionRelative
## 1 298.57557 2574023 1204246 -291.92704 -62.65759 W
## 2 109.16573 2574354 1204410 39.13835 101.90852 N
## 3 348.02859 2574052 1204080 -262.57843 -228.42169 W
## 4 49.28907 2574305 1204356 -10.44118 48.17047 N
## 5 266.77146 2574170 1204532 -144.92578 223.97216 N
## 6 858.50801 2573869 1203575 -446.42640 -733.30722 S
## closest_schreck tripIDnight DatetimeUTC hour day
## 1 WSS_2016_05 Ueli_110 2016-06-15 02:15:13 2 3Nachtviertel
## 2 WSS_2016_05 Ueli_114 2016-06-18 22:00:08 22 Abenddaemmerung
## 3 WSS_2016_05 Ueli_119 2016-06-24 02:15:08 2 3Nachtviertel
## 4 WSS_2016_05 Ueli_124 2016-06-28 22:15:11 22 1Nachtviertel
## 5 WSS_2016_05 Ueli_127 2016-07-01 22:15:07 22 1Nachtviertel
## 6 no_Schreck_witin_400m Ueli_133 2016-07-08 01:00:09 1 2Nachtviertel
## approachingRateRelativeS accelerationS sinousityS distanceS speedDiff3S
## 1 0.84857694 0.8251928 0.02771977 0.26834924 1.0000000
## 2 0.91306789 0.4437229 0.15740064 0.02366952 0.6761010
## 3 0.26017659 0.7787115 0.35117899 0.16542332 0.8457903
## 4 0.07122175 0.5537757 0.10448318 0.02006714 0.6912191
## 5 0.43267075 0.7000000 0.38061417 0.13690445 0.8276514
## 6 0.27509602 0.6354962 0.00000000 0.31355264 0.9759089
## approachingRateRelativeSlead scareEffect scareEffectShow scary
## 1 0.14715732 1.3819665 1.3819665 FALSE
## 2 0.01831587 0.9204379 0.9204379 FALSE
## 3 0.03056552 1.0773339 1.0773339 FALSE
## 4 0.09660522 1.0238393 1.0238393 FALSE
## 5 0.04343734 0.9666954 0.9666954 FALSE
## 6 0.01189577 1.2859567 1.2859567 FALSE
## return return_hrs highest scared_away
## 1 2016-06-18 21:30:08 91.24861 yes no
## 2 2016-06-24 01:30:11 123.00056 yes no
## 3 2016-06-28 21:30:06 115.49833 yes no
## 4 2016-07-01 21:30:40 72.00944 yes no
## 5 2016-07-07 22:30:08 144.99111 yes no
## 6 2016-07-11 21:45:18 95.25278 yes no ph<-ggplot(tt,aes(x=return_hrs))+geom_histogram(binwidth = 6)+theme_bw()+xlab("hours until animal returned within 400 of schreck [m]")
ph
### plot return hrs for scared versus non-scared:
sp<-ggplot(tt,aes(y=return_hrs,x=scared_away))+geom_boxplot()+theme_bw()+xlab("scared away")+ylab("hours until return")
sp 
#### scared away versus not
meanreturn1<-aggregate(return_hrs~scared_away, data=tt,mean)
nreturn1<-aggregate(return_hrs~scared_away, data=tt,length)
sdreturn1<-aggregate(return_hrs~scared_away, data=tt,sd)
## calculate standard error from sample size (n) and standard deviation (sd)
meanreturn1$SE<-sdreturn1$return_hrs/sqrt(sdreturn1$return_hrs)
### linear mixed model (accounting for repeated measures from the same individual)
m1<-lmer(return_hrs~scared_away +(1|TierName),data=tt)
summary(m1)## Linear mixed model fit by REML ['lmerMod']
## Formula: return_hrs ~ scared_away + (1 | TierName)
## Data: tt
##
## REML criterion at convergence: 735.7
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.3755 -0.3179 -0.2191 0.0967 5.0750
##
## Random effects:
## Groups Name Variance Std.Dev.
## TierName (Intercept) 9942 99.71
## Residual 5391 73.42
## Number of obs: 65, groups: TierName, 4
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 111.05 51.91 2.139
## scared_awayyes -26.72 20.17 -1.325
##
## Correlation of Fixed Effects:
## (Intr)
## scared_wyys -0.098drop1(m1,test="Chisq")## Single term deletions
##
## Model:
## return_hrs ~ scared_away + (1 | TierName)
## npar AIC LRT Pr(Chi)
## <none> 761.08
## scared_away 1 760.84 1.7602 0.1846confint(m1) ## 2.5 % 97.5 %
## .sig01 36.8245981 220.27793
## .sigma 61.6105647 88.14008
## (Intercept) -0.6285283 228.16184
## scared_awayyes -66.6196117 13.10484# assumptions of normal distribution of residuals violeted
# Main issue: only seven trips for scare equal TRUE
#### due to non-normality of residuals: non paratmetric test instead (however, not controlled for individual)
wilcox.test(return_hrs~scared_away,data=tt)##
## Wilcoxon rank sum test with continuity correction
##
## data: return_hrs by scared_away
## W = 490, p-value = 0.5745
## alternative hypothesis: true location shift is not equal to 0####
meanreturn<-aggregate(return_hrs~TierName, data=tt,mean)
nreturn<-aggregate(return_hrs~TierName, data=tt,length)
sdreturn<-aggregate(return_hrs~TierName, data=tt,sd)
## calculate standard error from sample size (n) and standard deviation (sd)
meanreturn$SE<-sdreturn$return_hrs/sqrt(sdreturn$return_hrs)
meanreturn$lwr<-meanreturn$return_hrs-meanreturn$SE
meanreturn$upr<-meanreturn$return_hrs+meanreturn$SE
pm<-ggplot(meanreturn,aes(y=return_hrs,x=TierName))+geom_point()
pm<-pm+geom_errorbar(data=data.frame(meanreturn),aes(ymin=lwr,ymax=upr),width=0.4)+theme_bw()+xlab("animal ID")+ylab("mean and SE of hours until return")
## boxplot
ptt<-ggplot(tt,aes(y=return_hrs,x=TierName))+geom_boxplot()+xlab("animal ID")+theme_bw()
ptt
## linear model
mp<-ggplot(tt,aes(y=return_hrs,x=scareEffect))+geom_point()+geom_smooth(method="lm")+theme_bw()+ylab("hours until return")
mp
m<-lmer(log(return_hrs)~scareEffect +(1|TierName),data=tt) ## TierName as random effect to control for repeated mesures per animal
qqnorm(resid(m)); qqline(resid(m))
summary(m) ## Linear mixed model fit by REML ['lmerMod']
## Formula: log(return_hrs) ~ scareEffect + (1 | TierName)
## Data: tt
##
## REML criterion at convergence: 116.7
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.7962 -0.5073 -0.1960 0.0249 4.8149
##
## Random effects:
## Groups Name Variance Std.Dev.
## TierName (Intercept) 0.7508 0.8665
## Residual 0.2915 0.5399
## Number of obs: 65, groups: TierName, 4
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 4.3731 0.5249 8.331
## scareEffect -0.2923 0.2326 -1.257
##
## Correlation of Fixed Effects:
## (Intr)
## scareEffect -0.531drop1(m,test="Chisq")## Single term deletions
##
## Model:
## log(return_hrs) ~ scareEffect + (1 | TierName)
## npar AIC LRT Pr(Chi)
## <none> 123.72
## scareEffect 1 123.35 1.6336 0.2012### correlation
cor.test(tt$scareEffect,tt$return_hrs, method="spearman") # spearman because of non-normality of residuals##
## Spearman's rank correlation rho
##
## data: tt$scareEffect and tt$return_hrs
## S = 49177, p-value = 0.5544
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
## rho
## -0.07467383